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A SUPERCONVERGENT HDG METHOD FOR THE INCOMPRESSIBLE 
NAVIER-STOKES EQUATIONS ON GENERAL POLYHEDRAL MESHES 

WEIFENG QIU AND KE SHI 


Abstract. We present a superconvergent hybridizable discontinuous Galerkin (HDG) method for the 
steady-state incompressible Navier-Stokes equations on general polyhedral meshes. For arbitrary conforming 
polyhedral mesh, we use polynomials of degree k-\-l, k to approximate the velocity, velocity gradient and 
pressure, respectively. In contrast, we only use polynomials of degree k to approximate the numerical trace 
of the velocity on the interfaces. Since the numerical trace of the velocity field is the only globally coupled 
unknown, this scheme allows a very efficient implementation of the method. For the stationary case, and 
under the usual smallness condition for the source term, we prove that the method is well defined and that 
the global L^-norm of the error in each of the above-mentioned variables and the discrete i7^-norm of the 
error in the velocity converge with the order of /c -|- 1 for k > 0. We also show that for A; > 1, the global 
L^-norm of the error in velocity converges with the order of k-\-2. From the point of view of degrees of free¬ 
dom of the globally coupled unknown: numerical trace, this method achieves optimal convergence for all the 
above-mentioned variables in L^-norm for k > 0, superconvergence for the velocity in the discrete i7^-norm 
without postprocessing for fc > 0, and superconvergence for the velocity in L^-norm without postprocessing 
for k > 1. 


1. Introduction 


In this paper, we consider a new hybridizable discontinuous Galerkin (HDG) method for the steady-state 
incompressible Navier-Stokes equations, which can be written as the following first order system: 


(l.la) 

L = Vm 

in H, 

(1.1b) 

—i/V ■ L - 1 - V • (u 0 m) -f Vp = / 

in H, 

(1.1c) 

V ■ M = 0 

in H, 

(l.ld) 

u —0 

on dQ 


(Lie) f P = 0 , 

Jn 

where the unknowns are the velocity u , the pressure p, and the gradient of the velocity L. v is the kinematic 
viscosity and / £ L^iVl) is the external body force. The domain C d = 2 , 3 is polygonal (d = 2 ) or 
polyhedral (d = 3 ). 

To define the HDG method, we adopt the notations and norms used in [6]. We consider conforming 
triangulation 7 )i of D made of shape-regular polyhedral elements which can be non-convex. We denote by £h 
the set of all faces F of all elements K G Th and set dTh ■= {dK : K G Th}- For scalar-valued functions cj) 
and '0, we write 


KeTh KeTh 

Here denotes the integral over the domain D C and {■,-)d denotes the integral over D C 

For vector-valued and matrix-valued functions, a similar notation is taken. For example, for vector¬ 
valued functions, we write (0),■0)7-^ := matrix-valued functions, we write (0,0)7-^ := 

T l<i,7<n i4'ij J '4’ij)Th ■ 
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Like all other HDG schemes, to define the HDG method for the problem, we introduce an additional 
unknown numerical trace which is the approximation of the velocity on the skeleton of the mesh. Namely, 
our HDG method seeks an approximation (Lh,Uh,ph:Uh) G Gh x Vh x Qh x to the exact solution 
(L| 7 y, u\eh) where the finite dimensional spaces are defined as: 

G^ := {G e L2(D) : G\k €^k{K), ViL € 

Vh := {v G L\n) : v\k G Pk+i{K), VX G %}, 

Qh ■■= {p & Ll{Vi) ■. p\K&Pk{K), MK &%}, 

Mh := {/x G L\£h) : v\f G Pk{F), WF G £h}, 

■= '■ v\da = 0 }. 

Here Pi{D) denotes the set of polynomials of total degree at most I > 0 defined on D, Pk{D) denotes the set 
of vector-valued functions whose d components lie in Pk{D), Pk{D) denotes the set of square matrix-valued 
functions whose d x d entries also lie in Pk{D), and Lo(^) = {P ^ L‘^{£l) : f^p = 0}. 

The method determines the approximate solution by requiring that it solves the following weak formula¬ 
tion: 

(l-2a) (L^, G)r^ + (uh, V • G)r^ - (uh, Gn)gri, = 0, 

(1.2b) (iyLh,Vv)rf, - {uh ® Uh,Vv)r^ - {ph,V ■ v)r^ - {yPhn-pun - {uh®Uh)n, v)gr^ 

-(^(V • Uh) Uh, ■y)rfc + {^{uh 0 (uh - Uh))n, v)oth = if, t>)Th, 
(l-2c) -{uh,Vq)THF [uh-n,q)dr^=Q, 

(1.2d) {vZhU-phU- {uh®Uh)n,pi)drH = 0, 

for all (G, V, q, fi) € Gh x Vh x Qh x M°. Here, 

(1.2e) {nLh-ph)n:=h'Lhn-phn-'^{nMUh-Uh)-Tciuh)iuh-Uh) on 8%, 

(1.2f) Tc{uh) := max(u/^ ■ n, 0) on each F G dTh- 

Here Um is the projection onto Mh- The goal of this paper is to consider the analytical aspects of the 
method including the rigorous proof of the uniqueness and existence of the solution of the above system 
and the error estimates for all unknowns. The computational aspects of the method will be discussed in a 
separate paper. 

Like other HDG schemes, the HDG method (11.21) uses the numerical trace of the primary variable Uh 
as the only globally-coupled variable. Our formulation is close to that of the HDG method in O [23], in 
which they have the same global degrees of freedom Uh- Nevertheless, there are three crucial differences 
which lead to special properties of our HDG method. Firstly, our method uses Pk+i{Th) to approximate 
the primary variable, which is the velocity, on each element while methods in [6l|23] uses Pk{Th) instead. 
Secondly, the diffusion part of the stabilization function ^{IlMUh — Uh) in (I1.2gI) is totally different from 
those used in [6l|23]. Finally, motivated by the work in [29] and |8], we insert two terms —(^(V • Uh) Uh, 
and {\{uh 0 [uh — Uh))n, v)gj-^ in (ll.2b|) . As a consequence of the above modifications, the new HDG 
method allows to use general polygonal mesh with optimal approximations for all unknowns. In addition, 
from the implementation point of view, like many other numerical methods for Navier-Stokes equation, we 
apply the classical Picard iteration to obtain the numerical solution. In each iteration we need to solve an 
Oseen equation. Due to our modification in (ll.2b|) . we can use the convection field obtained from previous 
step directly without the use of postprocessing. We notice that by adding these two terms in (|1.2b|) . the 
HDG method (11.21) is not locally conservative (see [5] for detailed explanation). On the other hand, since 
the HDG method (11.21) has high order accuracy for the approximations to all variables, lack of being locally 
conservative is only a minor issue. 

It is worth to mention that the HDG method using an enhanced space for the primary variable and the 
special stabilization function like ^{IlMUh — Uh) was first introduced by Lehrenfeld in Remark 1.2.4 for dif¬ 
fusion problem in [20] . He numerically showed that the methods provide optimal order of convergence for all 
unknowns without analysis. In [25] . we gave rigorous analysis for this approach for linear elasticity problems. 
Optimal order of convergence for all unknowns is obtained for both equations. In [24] . Oikawa analyzed a 
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HDG method for diffusion problem which uses the same polynomial spaces as in m, with a different choice 
of the numerical flux, he proved the optimality of the method for all unknowns. Since the polynomial order 
of the numerical trace of these HDG methods coiinniiii] is one less then that of the approximation space 
of the primary variable, from the point of view of degrees of freedom of the globally coupled unknown, they 
obtain superconvergence for the primary variable without postprocessing. In addition, all these methods 
work on general polyhedral meshes. However, the standard stability analysis for these methods can only 
provide the upper bound of 


{\\S/Uh\\'^^+h ^WUmUh — UhWgj-^) ^ , 

which can not control the standard discrete H^-norm of u^. We would like to emphasize that the control of 
the standard discrete H^-norm of Uh is essentially necessary in the proof of the HDG method (II.2|) having 
a unique solution. In [55], roughly speaking, we prove that 

(1-3) {\\VUh\\\+h ^\\Uh — UhWgj-^) ^ < C[WVUhW'^^ + h ^WllMUh — UhWgj-^) ^ . 

The inequality (11.311 enables us to control the standard discrete H^-norm of of the HDG method (|1.2I) 
such that we can show our method has a unique solution under the usual smallness condition for the source 
term. 

In literature, there are many existing mixed and DG methods designed for Navier-Stokes equations. See 
the classic mixed methods iniiKia, the stabilized methods proposed in [T2l[T6l[T9] and the DG methods 
[n (mill iniEHi EnisiiHiiii dni m] . An IP-like method and a compact discontinuous Galerkin (CDG) method 
were introduced in [22]. Recently, In [3], a mixed method is developed based on the pseudostress-velocity 
formulation. This method uses row-wise fc-th order RT space and Pk{Th) to approximate the pseudostress 
and velocity, respectively. The method provides the convergence of both variables in norm of order 
k + 1 — In|4], a modified method is introduced and it approximates the pressure directly with same 

convergence rate as in [3]. Later in a new mixed finite element method is introduced in which the stress 
is the primary variable. This method uses Pfe(7/j), row-wise fc-th order RT space and Pk{Th) to approximate 
the velocity gradient, stress and velocity, respectively. The convergence of the velocity gradient and velocity 
in L^-norm and the stress in iJ(div)-norm is of order fc. More recently, in 2015, Gockburn et al |5] gave an 
error analysis of the HDG method developed in [23] which is close to method in this paper. Our method may 
be criticized by the fact that with the modification of the scheme, we no longer have the local conservation 
of the momentum. Nevertheless, our approach has several advantages comparing with the one in [61123]. For 
instance, the analysis in [6] [23] is only valid for simplicical meshes and it needs a postprocessing procedure to 
obtain superconvergent approximation to the velocity. The method in [ida does not have superconvergent 
approximation to the velocity in the discrete iJ^-norm if fc = 0, while our method does even this the lowest 
order case. From the implementation point of view, in each iteration, the scheme in [idi] needs to solve 
a Oseen equation using a postprocessed convection field from the previous iteration while in our scheme we 
can use the velocity field directly from the previous iteration without any postprocessing. 

In this paper, we prove that the discrete iL^-norm of the error in the velocity, the L^-norm of the error 
in the velocity, the pressure and even in the velocity gradient converge with the order fc -|- 1 for any fc > 0; 
and that the velocity, for fc > I, converges with order fc -|- 2. Notice that as a built-in feature of HDG 
methods, see m, the degrees of freedom of the globally-coupled unknown comes from the numerical trace 
of the velocity on the mesh skeleton. From the point of view of the global degrees of freedom, the method 
provides optimal convergent approximations to the velocity, velocity gradient and pressure in L^-norm for 
fc > 0, superconvergent approximation to the velocity in the discrete iL^-norm without postprocessing for 
fc > I, and superconvergent approximation to the velocity in norm without postprocessing for fc > 1. In 
addition, the analysis of our method is valid for general polyhedral meshes. To the best of our knowledge, 
no other known finite element method for the Navier-Stokes equations has all of these properties. 

The rest of paper is organized as follows. In Section 2, we introduce our HDG method for the problem 
and present the main a priori error estimates. In Section 3, we present some preliminary inequalities and 
stability estimates. In Section 4, we prove the existence and uniqueness of the numerical solution. In Section 
5, we provide the detailed proof of the main results. 
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2. Main Results 


In this section, we present the main error estimates results. To state the main results, we need to introduce 
some notations. We use the standard definitions for the Sobolev spaces for a given domain D with 

norm 

\a\<e 


For vector- and matrix-valued functions (f> and $, we use ||</>||i,p,_D = ll^ilk.p.D, and ||$||i,p,D = 

j=i Moreover, whenp = 2 and £ < oo, we denote W^’^(L)) by H^{D) and || • ||f, 2 .D, by || • \\e,D- 

When I = 0,we denote W^’^{D) by L^{D) and the norm by || ■ ||lp(d): when I = £) and p = 2, we denote the 
L^{D) norm by || • Hu. 

We also introduce the following norms and seminorms: 


|||(t’,M)lllo,/t (ll^’llrfe + (II^A MllaTh + (^“M)llarh)) 

|||(’^,m)|||oo,/« := ||’^||L~(n) + IImIIl~(£^) 

Here || • ||ar^ := {Y^KeU H ’ ■ We also set 


yiv,fi) in H\Th)xL\Sh), 

y{v,fi) in H\rh)xL^{£h), 
V iv,n) in X L°°{£h)- 






where the average of v, f vj, is defined as follows: On an interior face F = dK~ n dK~^, we have f := 
+ v~)^ where denote the trace of v from the interior of and is the outward unit normal to 
K^. On a boundary face F C dK~ O dil, we formally set := v such that = v on dfl. We note that 
II • lli.ii is the standard discrete Ff^-seminorm. 

We are now ready to state our first main result on the existence and uniqueness of the numerical solution. 


Theorem 2.1 (Existence, uniqueness and stability). If ||/||a is small enough, the HDG method (121) has a 
unique solution (L, Uh,Ph, Uh) € Gh x Vh x Qh x Furthermore, the following stability bound is satisfied 

(2.1) |||(M,,,ti,,)|||i,?i < Ou-^II/IIq. 

for some constant C independent of v, the discretization parameters and the exact solution. 

Next we present the error estimates result for all unknowns. In order to have optimal error estimate 
for the velocity, we need some regularity assumption of the following dual problem. Consider the problem 


of seeking («/>, '0) such that 

(2.2a) $ - V</. = 0 in H, 

(2.2b) -i/V-$-V-(</>(8)M)-Vi/^-i(V</>)^u-f i(VM)^^ = 0 in H, 

(2.2c) V • </. = 0 in H, 

(2.2d) (p = 0 on d£l. 


Assume that the solution to the dual problem satisfies the following regularity estimate: 

(2.3) ||$||i,n + ||0||2,o + ||V'lli.o<a||0||o. 

Remark 2.2. If ||tt||ui(n) is small enough compared with the diffusion coefficient u, the dual problem (j2.2ll 
has a unique solution ((/),'0) S Hq{£1) x il^(H)/]R. In fact, when we use the standard energy argument, we 
need to have 

(2.4) ^l((V«)T0 - 0)n| < u||V0||^ 

to obtain energy estimate of </>. We notice that l|((Vit)(V0)^ti, (/))o| < C||M||ui(n)||^|lui(n). It is easy 
to see that (12.4|) holds if ||w||ui(a) is small enough compared with the diffusion coefficient u. This completes 
the proof of the above claim. If we further assume u £ W^’^(H) (~l i°°(H), then, the regularity assumption 
(1231) comes from a standard regularity estimate m for the Stokes equations. 
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Now we are ready to present our second and main result: 

Theorem 2.3. If ||/||n is small enough, then we have 

||L -L/illn + \\u-Uh\\Q + \\u-Uh\\i,h + \\p-Ph\\Q < Ch^+^, 

Here the constantC depends on ||M||fc+ 2 ^n, ||p||fc+i_n,u and k. In addition, if the regularity assump¬ 
tion (lOl holds and u G then for k >\ we have 

\\u - uiilln < Coh^^^- 

HereCn depends on ||u||i«.(o), ||ri||fc+2,n, ||p||fc+Ln, u and k,Cr- 

3. Preliminary estimates 

In this section, we present some preliminary inequalities for the proof of our main results. First, we would 
like to recall an important inequality which was introduced in |26j . Here we write it in a slightly general 
way. Though our results in this section and the following ones are valid for conforming meshes with shape 
regular assumption, we assume the meshes are quasi-uniform for the sake of simplicity. 

Lemma 3.1. For any given function (L, v, /x) G Gh x Vh x satisfying p.2ap . then we have 

IIK'*^, < CffDodlLlIn + h~^\\nMV- fiWdTh)- 

For the proof of the above result, we refer the Lemma 3.2 in [26]. In addition, we also need the following 
basic inequalities: 

Lemma 3.2. For 1 < q < oo {d = 2), 1 < g < 4 (c? = 3), there exist positive constant Cq such that 
(3.1a) ||t»||L.(n) yvGV{h), 

(3.1b) ||t,|U,(n) < C',|||(«,m)IIIi.i., V (ti,/x) e Vih) X i^. 

Here V{h) := iJg(H) + Vh- In addition, we have a trace inequality: 

(3.1c) ||HlL4(ar.) < Ch-^\v\\,^h < Ch-^{v,^^)l,^h. V (ti,/x) e V{h) x 

The proofs of (I3.1all - (l3.1c|) are provided in Proposition A.2 in |6], Proposition 4.5 and (7.7) in [18]. To 
simplify our notations, we group all the nonlinear terms in our formulation as the following operator: 

Definition 3.3. For any {w, w), {u,u), {v,v) G H^{Th) x L^{£h), we define the operator: 

0{{w,w)-, {u,u), {v,v)) := - {u®w,Vv)r,, “ - w)n,v)Qry, 

+ {Tc{w){u-u),v-v)orh + {{u®w)n,v-v)QrK- 

The above operator plays a crucial rule in the analysis. It has the following coercive property: 
Proposition 3.4. For any {w, w), {u,u) G H^{Th) x L^{Eh), ifu\dvi = 0; then we have 
0{{w,w)-,{u,u),{u,u)) = {{tc{w) - ^w- n){u-u),u-u)9rH > 0. 

Proof. By Definition 13.31 we have 

0((w, w); {u,u),{u,u)) := - {u(^w,Vu)r,^ “ w)n, u)arh 

+ {Tc{w){u-u),u-u)oTh + {{u(^w)n,u-u)orh- 

Applying integration by parts for the first term, we have 

(m® w,'s/u)rt, = -(V • (m® w), u)r,, -f ((m® w)n, u)drh 

= -((V • w)u,u)r^ - {u®w,Vu)r„ -f {{u®w)n, ujdTh- 

-(tt® w,Vu)r,, - (i(V • w)u,u)r,, + {[]^u®w)n, u)or„ = 0. 


This implies that 
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Inserting above identity into 0{{w, w)-, {u,u), {u,u)), we have 

0{{w,w);{u,u),{u,u)) = {Tc{w){u-u),u-u)on + {{u<S:w)Ti,u-u)arH - {^{u0w)n,u)or^ 

= ■ n){u - u), u - u)or^ - {^{w ■ n)u,u)arh 

= {(Tciw) - ^w- n){u-u),u-u)oTH > 0- 

The last step is due to the fact that u is single valued on £h and M|ao = 0. □ 

Next, we present a continuity result for the nonlinear operator O that we will use throughout the analysis. 
We first define the following space: 

Hl{^) := {(u;,S;) e Hl{£l) x L\£h)\w\E, = w}, 

The above space is the graph space of the trace mapping from onto L^{£h)- We are ready to state 

the following result: 

Lemma 3.5. There is a positive constant Co such that 
(3.2) 

\0{{wi, Wi)] (u,u), {v,v)) - 0((W2, W 2 ); iu,u), (u, u))| < Coiiwi, Wi) - {W 2 , ti> 2 )|||i,/» |||(m, IK?;,9)111,/,, 

for all {wi, Wi), {w 2 , W 2 ), (u, u) G ffo(^) + (^h x -^) and any {v,v) G Vh x 
Proof. Setting 6^ ■= Wi — W 2 , <5® := Wi — W 2 , by the definition of the operator O, we have 
0{{Wi,Wi)]{u,u),{v,v)) - 0 {{W 2 ,W 2 )] {u,u), {v,v)) = 

- {u(S)Sy„Vv)rh - (^(v • Sy,)u, v)rf, + (i??® (<5^ - S:s,)n,v)oTH. 

+ {(.Tc{wi)-Tc{w2)){u-u),v-v)on + {{u^S^)n,v-v)QrH, 

applying intergration by parts in the first term, rearranging the terms, we have 

= v)r^ + {v0 6y„\7u)r^ + (i?x® {6n, - 5^)n,v)oTh 

+ {{tc(.wi) - Tc{w 2 )){u-u),v-v)or^ + {{u <^1 -u®5yf)n, v)arh 

= Ti + Ta + T 3 + Ti + ^5. 

Next we estimate each Ti. 

For Ti, we apply the generalized Holder’s inequality, 

Ti < ||V-(5™||rJ|w|U4(0)||?;|K4(Q) < C|||(l5^„,(5s,)|||l.h |||(m, w)|||i.;i |||(v, v)|||i.h, 

the second inequality is due to (I3.1bl) . T 2 can be bounded in a similar way. 

For Ta, we apply the generalized Holder’s inequality, 

Tz < h~^\\5n,-5^\\L^{dTh) h'^WAL^tdTh) < C'lll(<^«>:<^S,)|||l,/i lll(“: 

the second inequality is by (I3.1c|) . 

For T 4 , by the generalized Holder’s inequality we have: 

Ti<C\\Tc{wi)-Tc(w2)\\L'^(^orH) \\'^~M\L^(dTH) ll«- 

by the fact that the function max(??; • n, 0) is Lipschitz, 

< C\\5is,\\L4.(^arh) \\'*^-'^\\L'^{dTh) 

< + \\5w\\L*(dTh)) lll(“’“)llliA ll’^-^llL'‘(arh)) 
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Notice here if (Sy,, (5s,) G + (Vh x then (5^ — i5sj G VhWh- Hence, we can apply inverse inequality 

on ||(5«, - (5s,||i4(ar;,), llu-911^4(97-^) to have 


<Ch^{h 4 ||(5«, - i5s,||l 2(97;) + ||5^||i4(97-^)) 

< C'/i3|||((5„,5si)|||i,,, |||(M,tt)|||i,h IK?;, ?))|||i,h, 


tt, u 


\\l,hh 4 11 ^- 9 ||i 2 ( 97 -,^) 


by p.lcl) and the fact that d = 2, 3. 

Finally, for T5 we first break it into two terms: 


Ts = -{{uiSi {6y,-6^))n,v)or^ - {{{u-u)(^Sy,)n,v)Qr^. 

For the first term, by the generalized Holder’s inequality, we have: 

(( m (8 iSy,-dji,))n,v)QrH < h^imiL^dn) h~i\\6y,-Si,\\L2(0TH) 

< h^{\\u- u||l 4(97^) + ||tt||L4(97-,^)) h~2\\Sy, - (5^11^2(97-^) WAL^^idTh) 

<C|||((5.,^s,)ll|i,. I|IK«)|||i,k III(«,«)|||i,k. 

In the last step we used the same argument as in for the term ||m— m|| 2,4(97-,^) and (|3.1cF The second 
term can be estimated in a similar way. We complete the proof by combining the estimates of T^. □ 


4. Uniqueness and existence of the numerical solution. 

We will apply the Picard fixed point theorem to show the existence and uniqueness of the solution of 
(O. To this end, we begin by rewriting the method into a more compact and appropriate form for the 
proof. If we add (ll.2aD - (II.2dl) . the method can be written as: Find (L/i, Uh,Ph,Uh) £ Gh x Vh xQhxM^ 
such that 

(4.1) S{{Lh, Uh,ph, Uh), ((G, V, q, fi))) + 0{{uh, Uh); (uh^Uh), (v, /x)) = (/, v)r^, 
for all (G, V, q, /x) G Gii x Vh x Qh x M°. Here the bilinear form 5(-, •) is defined as: 

S{{Lh,Uh,Ph,Uh), (G, i;,g,/x)) := {Lh,G)u + (ix^, V • G)u - {uh,Gn)oTh - (u V • pL,,)r,, + {n,vLhn)oT^ 

- {uh,Vq)rh + {uh ■ n,q)Qrh + {v,Vph)Th “ ’ 'n,Ph)dTh 

+ {-^{IlMUh - Uh), V - p)dTh- 

We also define a mapping T as follows: for any {w, w) G H^{Th) x L^{£h), {uh, Uh) = J-{w, w) £ VhX Mh 
is part of the solution {Lh, Uh,Ph,Uh) £ Gh x Vh x Qh x of 

(4.2) S{{Lh, Uh,Ph,Uh), {G, v,q, p,)) + Oi{w,w); {uh,Uh), {v, p,)) = {f,v)r^, 

for all (G, V, q, p) £ Gh x Vh x Qh x M^. It is worth to mention that when w £ H{div, fl), V ■ w = 0 and 
w = w\£f^, the above system is a HDG scheme for the Oseen equation. Clearly, (uh, Uh) is a solution of (11.21) 
if and only if it is a fixed point of the mapping J-. Next we present a stability result for the above scheme. 

Lemma 4.1. If (Lh, Uh,Ph,Uh) £ Gh x Vh x Qh x is a solution of (14.2F then there exists a constant 
C that solely depends on the constants Chdg o,nd C 2 in Lemma \S. 1[ Id.Sl such that 

Proof. Taking (G, v, q, p) = (uL/j, Uh,Ph, Uh) in (14.21) . after some algebraic manipulation, we have a simplified 
equation: 

^l|L/i||n + {^{nMUh - Uh), Uh - Uh)dTh + 0((w, w); (uh,Uh), {uh,Uh)) = (/, Uh)n. 

Therefore, by Lemma IXTl Young’s inequality and Proposition 13.41 we have 

< 2CHQQu(||L;i||Q + — \\IlMUh — Uh\\gj-I^) 

< 2C'hdg ('^IlL/illn + {^(nnUh - Uh), Uh - Uh)aTh + 0((w, w)-, {uh,Uh), (uh,Uh))'^ 

= 2Ghj3q(/, ti?i)n < 2C'HDGll^/^llf2|l/llf2 — 2G2 C'hdq|||(iX/^, ti;i)|||i7^||/||n. 


The last step is by Lemma [3.21 with q = 2. This completes the proof with C = 2C'2C'y]9Q 


□ 
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Inspired by the above stability result, we define a subspace of Vh x 

'■= & Vh X < C'2 C'hdg^ ^||/||n}- 

We are now ready to give the proof of the existence and uniqueness result for the HDG scheme ira/gii). 

Proof, of Theorem 12.11 

Clearly, T maps Vh x into ICh hence it maps ICh into itself. In order to show the existence and 
uniqueness of the solution of (fT^/dO). it suffices to show that is a contraction on the subspace ICh- 
To this end, let {wi,wi),{w 2 ,w 2 ) £ Kh and (L^, are the solutions of the problem (I4.2|) with 

{w,w) = {wi,Wi), {i = 1,2). So we have := and := J^{w 2 ,W 2 ). If we set 

(5l := Li — L 2 , 5u '■= Ui — U 2 , Sp := pi —p 2 and 6% := U 1 — U 2 , due to the linearity of the operator 5, we have 

>5(((5l, Su, Sp, S^), (G, V, q, n)) + 0{{wi, wi)-,{ui,ui), {v, n)) - Oi{w 2 , m); (m 2 , M 2 ), (m, m)) = 0, 

for all (G, V, q, fi) £ Gh x Vh x QhX M°. Taking (G, v, q, pi) = Sp, 5%) into the above identity, after 

some algebraic manipulations, we obtain 

^^II^LlIn + {^{HmSu- S^),6u - SigfdU = -0{{Wi,Wi); (Mi,Mi), (5„,(5b)) + 0{{W2,W2)-, (^ 2 ,^ 2 ), {Su,Sn))- 
Or 

^\\Sh\\n + {^{HmSu - Su), Su - Si}dTh+0{{w2, 102 ); (<5«,^b), (^«, ^s)) = 

0{{W2, M>2); (mi. Ml), (Sn,S^)) - 0((wi, Wi); (mi. Mi), ((5„,5 b)). 

By Lemma [3TT] Young’s inequality and Proposition 13.41 we have 

^M)||ln/i < 2GHQQ(z^||5L||n + — S%\\qj-^) 

^ + {-^{IImSu — ^s), Su — SiffoTh + 0[[w2, W 2 )', (Su, S^), (Su, ^s))^ 

= 2 Ghdg ((^((‘u>2, W 2 ); (uuui), (Su,S^)) - 0((wi, wi); (ui,Ui), (Su,S^))^ 

< 2GoC'hdgIII(^i - '*^ 2 , m - W 2 \\\i,h |||(mi, Mi)|||i,,i |||(^u,<5s)|||i,/t by LemmallSl 

< 2C'ciG2C'HDG^^~^ll/lln IIK^^I - 1«2,M>1 - W 2 \fii,h |||(^u,(^s)|||i,/i by Lemma|4Tl 

Therefore, we have shown that 

lll(^«,<52)|||l,;i < 2CciG2C'HDG'^”^ll/lln III(m^1 - W2,Wi - W2\\\i^h- 
Obviously, the above bound implies that is a contraction on ICh equipped with ||| • |||i/j provided 

II ill / 

11 ,71 1 O — r)/^ /^4 

20 c)G2Gjjj3q 

By the fixed point theorem, there is a unique fixed point (uh,Uh) £ fCh of the mapping T. It is also the 
unique solution of the system (O. This completes the proof. □ 

5. Proof of the error estimates 

In this section, we provide the detailed proof of the main error estimates for all unknowns. We proceed 
in several steps. 

Step 1: Error equations. We begin by introducing the error equations that we are going to use in the 
analysis. For convention, we introduce the following notations for the errors: 

Cl := 77gL — L/i, e„ :=/7y m — m/i, Cp := Uqp — ph, := 77mm- m/j, 

(5l := L - Z7gL, (5„:=M-/7yM, Sp := p - Uqp, (5s:=m-77mm. 
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Here IIo, Uy, Hq, Um are the projections onto Gh, Vh,Qh, respectively. In addition to the basic 
inequalities listed in Lemma 13.21 we will frequently use the following basic inequalities as well: 


(5.1a) \\q\\F <Chj^^\\q\\K, 

(5.1b) \\D^{q - niq)\\K < 

(5.1c) \\q-mq\\F<Ch'^^\\q\\i+yK, 

(5.1d) \\q-nMq\\F<Ch'^+^\\q\\k+l.K 


for all q G Pi{K), (/ > 0), 
for all q G 0 < m < Z, 

for all q G 0 < m < I, 

for all q G 


Here 77; denotes the projection onto Pi{K), F denotes any face of K. In addition, we have the following 
estimate for the projections under the triple norm ||| • 


Proposition 5.1. For any uG we have 

(5.2) IIKTIyM, 77 mm)|||i,;i < CllttHi^n. 


Proof. By the definition of the norm ||| • we have 

|||(77yM, 77 mm)|||i,;i < 2(||V77yrt||7-^ + h~^\\IIvu — IlMu\\dj-h)- 
We are going to bound each of the above terms by ||M||i,n- For the first term, we have 

\\vnvu\\r^ = ||V(77yM- u)\\r^ < Ch~^\\nvu- u\\o,Th^ 


here u denotes the average of u within each element K G Ft, the inequality is by the inverse inequality of 
the polynomial spaces. 


by the Poincare inequality for each K G Th- 

For the second term, applying a triangle inequality we have 

h~^\\nvU- nMuWgr^ < h~^{\\u— UynWorh + l|w- nMu\\gr^) < 2h~^\\u- IlkuWgr^, 

here Ilk denotes the projection onto Pk{K) for each K G Th, 

< C'(||V(m— nkUjWq-,, + IlkuW-f-ff) by the trace inequality, 

the last step is by a similar argument as for the first term and (ED. This completes the proof. □ 

It is not hard to verify that the exact solution (u, L,p, rt|f:^) satisfies the following equation: 

5((L,it,p,«),((G, ■y,g,p))) +C>((tt,M);(it,M),(i;,p)) = + (ps, v- y)gr,^, 

for all (G, V, q, fi) G Gh x Vh x Qh x M^. Subtracting (14.IL we have 
5((L, u,p, u), ((G, V, q, p))) - <S((L;i, Uh,Ph, Uh), ((G, v, q, p))) 

+ 0{{U,U)- {U,U), (V,p)) - 0{{Uh,Uh); {Uh,Uh), (?I, p)) = (pB, V- p)gr,^, 


by the linearility of the first operator S, we have 

5((eL,e„, ep,eB), ((G, v, q, p))) + Oi{u,u); {u, u), {v, p)) - 0{{uh,Uh)\iuh,Uh), {v, p)) = 

^(('^L 5 ^b) j ((G, V, q^ A*))) “t” h’}dFh ' 

Finally, by the definition of the operator S and the orthogonality property of the projections we have 
the error equation: 

5((eL, e„, Cp, cb), ((G, v, q, fi))) + 0{{u, u); {u, it), {v, fi))-0{{uh, Uh)] {uh, Uh), {v, p)) = 

^ ^ {FShU-SpTi- -IImSu, V- p)grh, 


for all (G, V, q, fi) G Gh x Vh x Qh x iW°. 
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Step 2: Estimates for cl, e„. We first apply an energy argument to bound the errors en, e„, which 
can be stated as follows: 

Lemma 5.2. If the exact solution m, L,p is smooth enough, and ||M||i,n is small enough, we have 

||eu||n + |||(eu, e^)|||i,?i < CHDaiW^'LWn + h~^ Wnm^u - eigWdTh) < . 

Here the constant C depends on ||ti||fc+ 2 ,n, HitHoo.n, ||p||fc+i,n, and k but independent of h. 

Proof. Taking (G, v, q, /x) = (bl, Cu, Cp, e^) in the error equation (I5.3L the resulting equation can be simplified 
as 


or 


^"lleLlIn + ’^WnMeu 


^uWdTh + ■“)’ eu))-0{{uh, Uh); {uh,Uh), (e„, e^)) = 

^ dp 71 — Hm du , Cu C'^dTh 5 


i^lleLlIn + TWnMCu- CuWlrh =^ii'>J-h,Uh); {uh,Uh), (e„,ea)) - 0{{u, u); {u,u), (e„,es)) 

(5.4) * 

“f {^dl^ 71 dpTl "J^Hj^^du , Cu Cu}d'Th ■ 

Let us first estimate the last term of the above equation. Applying the Cauchy-Schwarz inequality, we have 
{vdi,7l — dpTl — -77m<5u , Cu — CisfdTh < (^ll<^L||9rh + ll'^plloTfe + vh~^\\nMdu\\dT,}j \\Cu — CuWdTh 

< (f^WdhWdTh + ||<5p||arfe + vh~^\\du\\aTh')^~^\\eu-CuWon 

< (^(II^lIIh + ll^plln + Il|(e„, eu)\\\i,h by (IS.laD . 

< C'/i''+^(||L||fc+i,n + ||p|U+i,n + ||M|U+ 2 ,n)|||(e„,es)|||i,?i by (IS.lbl) . 

The main effort in the analysis is to have an optimal estimate for the nonlinear terms, we rewrite these two 
terms into four parts: 


(^/l7 '^h)-i {c-ui Cu^^ ^((^7 7 (^7 7 {Cu, Cu')') 

— 0((t/./i7 Ufi), (tX/i, Ufi), {Cu, Cu)) l^(('li/i7 '^h) ^ {dlYTl, HjyU), {Cu, Cu)) . P± 

+ 0 {{uh,Uh)-, {nvu,nMu), {eu,eu)) - 0{{nvu, Hmu)] {nvu,nMu), (e„, e^)) : T 2 

+ OdHyu, Hmu); {Hyu, Hmu), (e„, e^)) - 0 {{u, u); {Hyu, Hmu), (e„, e^)) : T 3 

+ 0 {{u,u)-{nyu,nMu),{eu,eia)) - C>((m, m); (m, it), (e„, e^)) : T 4 


Notice that the operator O is linear with respect to the last two components, so we have 

— ^((l1/i7 '^h): (.Cu: Cu), (.Cu, Cu)) ^ 0 , 

by Proposition 13.41 For T 2 , we apply the estimate (13.21) in Lemma [?751 we have 

T2 < CoWliHyu, nMu)\\\i.h\\\(eu, Cu)\\\yh < || it|| ||| (e„, 6 ^) ||| ^ /^. 

For T 3 ,T 4 , if we directly apply the continuity property of the operator O (13.2|) . we will only obtain 
suboptimal order of convergence. We can recover optimal convergence rate by a refined argument for each 
term. We start with T 4 , by the linearity of O for the last two components, we have 


T 4 — ^(('^, '^), (du, du), (Cu, Cu)) 

— (du G It, ^Cu)'j\ ( 7 'c('^)(du du), Cu cif))((du G it)it, Cu c^Q-y,^ .— T41 T ^42 T ^43. 

We now apply the generalized Holder’s inequality to bound these three terms as follows: 

T 41 < ^ ||it||oo.A'||<5ii||A||Ve„||_fs- < Ch^^'^\\u\\ oo,0 II 1l|| A:+2,0 III (^U7 ^u) III l,/i 7 

KdTh 


^42 ^ ^ ^ II ^ ' "^11 oo,dKh 2 II du du II h ^ II Cu Cu || ^ f^ll^' ^l|ooA/i^ ll^ll^+l,S4|||(eu7 ^u) ||| l,ti, 

KGTh 

T 43 < C\\u- n||oo.£fc^''’^^||it|U+i,n|||(eu 7 e^)|||i,h. 
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For Fa, by the definition of O, we have 

^3 ^ ® ^^u)'Th ( 2 ^7 ^u^Th ^ ^ ) ^7 cITii 

+ {{rcillMu) - Tc{u)){nvu- nMu),en - eu)drh - {{IImu^ Su)n,eu - eu)dTh- 
Among the five terms in the above expression, the third term needs some special treatment in order to obtain 
optimal convergence rates. For the others, we can bound them in a similar way as for T 4 i,T 42 ,T 43 . For the 
sake of simplicity, here we show how to bound the last term and then focus on the third term. By applying 
the generalized Holder’s inequality on the last term, we have 

{{IImu® 5%)n,eu —eu)dTh ^ ^ \\nMu\\oo,dKh'^\\5u\\dKh~^\\eu — eis\\oK 

K&Th 

< II II m|| fc+i.n III (e„, e^) ||| 

In the last step we applied the inequality: 

IIHm^II oo, dK < C'll^tll oOjdK • 

This result can be obtained by a simple scaling argument. Finally, let us focus on the third term. We rewrite 
the term as follows, 

^ ^u)dTh ~ ^ 2 ® i^u ^w))^7 ^u)d7h ^ 2 ^ d'Th 

- {^{nvU®S^)n,eiB)dTh 

= T 31 + T 32 + r33. 

For r 3 i, we apply the generalized Holder’s inequality, 

T 31 < ^2 ll-^E^*||oo.aif/i^ ||<5„ — i5s||aR-/i“2 ||e„ — esllaif < Cft.^'''^||ii||oo,n||u||fc+i,n|||(e„, e^)|||i,/i. 

KdTh 

For T 32 , we have 

T 32 < ^ IlHyMlIocaif^”^ll^ullaiflle^llaif < C/i^''‘^||tt||oo,0111*11^+2,n(^^||es||arfi)- 
K€Th 

For T 33 , inserting a zero term {^{u^ Sig,)n, 6 ^) 97 ^ into T 33 we obtain: 

^33 = {7:{^u® 5u)n,eu)dTh < '^2 ll'^slloo.aA^”^ ||<^u||aR'(/i= ||es||aA) < ||w||oo,n||M||fc+2,n(/i2 lle^Horh)- 

KeTh 

The last step is to show that 

llsBlIaTl — C'll(e«, e-^) III. 

To this end, we apply a triangle inequality and (l3.1cF 

/i^lleslloTh < /»^l|e«||arfe + ||e« - esHor^ < C'(l|e«||a + /i|||(e„, ea)|||i,;i) < CIKen, e^)|||i,h. 

This completes the estimate for T 3 . Finally, if we combine the estimates for Ti — T 4 and Lemma I3.1I we 
obtain the result stated in the Lemma [521 D 

Step 3: Estimates for Cp. Next we present the optimal error estimate for Cp. As usual, we bound the 
pressure error via a inf-sup argument. It is well-known [5] that the following inf-sup condition holds for a 
polygonal domain Q: 

(5.5) sup for all g G Lg(H). 

™effJ(a)\{o} ll^lli.n 

Here k > 0 is independent of w,p. We can bound Cp using the above result. 

Lemma 5.3. Under the same assumption as in Lemma \5.2l we have 

\\ep\\a<Ch^+^. 

Here the constant C depends on ||ri||fc+ 2 ,n, ||r*||oo,n 7 ||p||fc+i,n 7 u, fc and k but independent of h. 
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Proof. Since tp € LQ{n), by (15.51) we have 


(5.6) 


llcplln < 


1 

— sup 

« w&HUljMO} 


(V ■ w, ep)n 

IIHIi.o 


Now let us work on the numerator. Applying integration by parts and the orthogonality property of the 
projections Uy, Um, we can rewrite it as follows: 


(V • w,ep)n = (ep, V • nvw)r^ + {{w- Flyw) ■ n,ep)oTh. = (cp, V • nywjr^ - {{Ilyw - Umw) ■ n, ep)ar^- 
If we take (G, v, q, /i) = (0, Ilyw, 0, Umw) in the error equation (15.31) . we obtain: 

(cp, V • nvw)rh - {{IIvw- Umw) ■ n,ep)dTh 

= VI7yiy)r^ — {veyn — -{IlMeu — es), Ilyw — nMw)dTh 
+ {v5i,n— Spn— —nMSu,nvw— nMw)gj-^ 

+ u)] (u, u), {IlyW, IImw)) - 0 {{Uh,Uh); iUh,Uh), {nyW^UMw))^ 

:= Ti + T2 + Ta + Ty 

Next we show that 

T^+T2 + T^ + Ti<Ch^+^\\w\\yn- 

For Ti we have 

Ti < y||eL||n||VI7yi(;||rh < Ch^~^^\\w\\i^n, 

by Lemma 15.21 and (EU. For T 2 , by the Cauchy-Schwarz inequality, we have 

T 2 < \\eh\\dTh + h~^\\IlMeu - eu\\drh)ih~^\\nvw - UmwIIoTh) 

<Ch^+^\\w\\yn. 


Ta can be estimated similarly as T 2 . Finally, for the last term, we split into three terms as: 

T4 = 0{{u, u)-{u, u), {Ilyw, IImw)) - 0{{u, u); {Uyu, IImu), {Uyw, IImw)) : r4i 

+ 0 ((m, m); (IlyM, Umu), {Ilyw, IImw)) - 0 {{uh, Uh)\{IIvu, Umu), {Uyw, IImw)) : T42 

+ 0{{Uh,Uh)\ {HyU, Hmu), {HyW, Hmw)) - 0{{Uh,Uh); {Uh,Uh), {nyW, Hmw)) : T 43 . 

We bound T 4 i separately. 

T 41 = 0 ((u, u); {Su, 6 %), {Ilyw, IImw)) 

= {Su <S> u, nyw)r,^ + {tc{u){Su - Su), IlyW - nMw)oTh + {{^u u)n, IlyW - nMw)gr^ 

< ||M||oo.n||^u||n||^Iyi(i||n + ||u||oo,n h^{\\Su — Su\\aTh + ll^ulloTfe) h~2\\nyw — IlMwllgr^ 

< Ch'"~^'^\\u\\^^Q\\u\\k+l,Q\\\{nyW,IlMw)\ji^h < G/l'"+^||w||oo,o||u)||l,n, 

the last step is by the approximation properties of the projections (|5.1ap , (IS.lbl) . 

For T 4 a, due to the linearity of O on the last two components and (13.2L we have 

T 43 = 0{{Uh,Uh);{eu,e^),{nyW, IImw)) < C\j{Uh,Uh)\h,h IIKe^, et;)|||y;i \j{nyw,nMw)\\\i^h 

< C'(|||(I7yM,I7MM)|||i,h + |||(e„,es)|||i,/i) |||(e„, e5)|||i,/i |||(77yu;, I7Mw)|||i,h 

< Ch>^+^\\u\\yn\\\{nyW,nMW)iyh < Ch’^+^Myn, 

by (El]) and Lemma ED 

For T 42 , if we directly apply EU, we will only obtain suboptimal convergence rates. Alternatively, we 
need a refined analysis for this term. First, we let Eu '■= u — Uh = e„ + and := u—u= eu + Su- Next, 
by the definition of O, we can write T 42 as 

T 42 = —{IlyU® Eu,VIIyw)j-i^ “ (^(^ ■ Eu)IIv'fJ', IIvw)j-h, + {^u “ Eu))n, 

+ {{tc{u) - Tc{Uh)){IIyU- nMu),nyW- nMw)gr,, - {{IImU® Eu)n, IlyW - IImw) gr,, 

= Si + --- + S 3 . 
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Notice that by Lemma 15.21 and (|5.1bl) we have 

(5.7a) ||L;„||i,^ < ||e„||i.,, + ||,5„||n^ < Ch>^+^\\u\\k+ 2 M, 

(5-7b) \\Eu\\dTh < + lle^llarh + l|e« - efiHar^ < , 

by (jS.ldl) . Lemma [5^ Now we bound each of Si. By the generalized Holder’s inequality, we have 
Si < ||77yu||oo,n||7?«||n||VTIvtallo < C'/i^''‘^||ii||oo,n||t<i||i,n. 

By a similar argument, we can bound S 2 as: 

For 54 , we apply the generalized Holder’s inequality to have 

Si < || tc ( m ) — Tc{uh)\\dTh\\nvu — nMu\\oo,dTh h~^\\nvw — IlMwWdTh 

<Ch^\\u— UhWdTh II w||oo,n h~^ Wllyw — UMwWaTh 
< C'h'=+^||M||oo,n||i9||i,o, 

by the estimate (I5.7b() and (15. ip . 

By a similar argument we can bound S '5 as 

S's < C'h''+^||M||oo,a||i(i||i,a. 

For the last term S' 3 , if we apply similar estimate as the others, we will only obtain suboptimal order 
convergence rates. Therefore, we need a refined estimate for this term. We rewrite S '3 as follows: 

(e« - e^))n, nvw)orh + 5u)n, nvw)Qr^ - {^{nvu(^ 6%)n, nvw)gr^ 

< ||77vu;||arfcl|77i/it||oo,9rh^~^(l|eu - e^llaTh + ll<^u||9rfe) - {-^{nvu® Sis.)n, nvw)drh 

< - ^{{nvU'»Su)n,nvw)oTh, 

by (j5.1cp . Lemma I 5 . 2 I For the last term, we further split it into two terms as: 

{{IIvu® 5u)n, nv'w)drh = {{IIvuiS: 6u)n, w)oTh - {{Uvu® 5u)n,w - nvw)orh 

= -{{{u-IIvu) ® 6is,)n,w)drh. - {{IIvu® 5^)n,w-nvw)drh^ 

where the last step is obtained by inserting a zero term {{u® 6%)n, w)g-j-f^ = {{u® 6%)n, w)gQ = 0 . 

< ^^ll<^>x||L'‘(arh) ^^ll^lU^(arfc) ll<^B||9rfe + ||77uM||oo,arfe||<^ji||9rfe||^«- nvw\\gr^ 

< C'/i''||i5„||i,;i||w||i,n + IIm|| 00 , nIIi«||i,n 

<C/i"+i(||ti|| 2 .Q + ||n|U.o)||Hluo, 

where in the last step we used the inequalities (I3.1cl) . (I5.1bl) . (I5.1dl) . The proof is complete if we combine all 
the above estimates. □ 

Step 4: Optimal estimate for e„. Notice that Lemmaprovides an optimal estimate for cl but only 
suboptimal estimate for e„. This is due to the fact that we use a Pfc+i polynomial space for the unknown u. 
To obtain an optimal convergence estimate for Cu we will use the adjoint problem (I2.2p to apply a duality 
argument. We begin by the following identity for the error eu- 
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Lemma 5.4. Let be the solution of the dual problem (12.21) with the source term 9 = e„, then we have 

II^ulln ^ d7h 

- {-r{nM&u — es), nvf) — nM<i))dTh 

h 

+ {v5i,n — SpU — —IImSu, Llycf — IlM^^dTh 

- ((e«, V • (0 (8> u))r^ + 0{{u, u); (e„, e^), {Hvcf), HmcI)))^ 

- 0((u, u); {Su, Su), {nv4>, Lluft))) 

+ {pifu, ; i.'^h ; '^h) 1 i^cj ); ^ 0 )) 7 "^/i) 7 7 ^/i); (^0 7 

+ (</>») -0{{u,u)\{uh,Uh),{(l>, 0)) - {eu,Y)r^ 

+ (z/CL, yS^)r^ - (cp, V • (5<^)rfe 
:= Ti + • • • + Tg. 

i/ere (K := and 

5cp\= (f) - nv<t>, 5^ ■='ip - IIqiI;, := (f> - nM4>- 
Proof. By the adjoint problem (I2.2al) - (I2.2cl) we have 

l|e«||a = - Heu, V • $)r^ - (e„, V • (0 0 - (e„, Vi/z)r^ - (e„, y)r^ 

- (i/eL, ^)u + '^(f))Th 

- (ep, V • 0)rfe 

rearranging the terms, we have 

= - z/(e„, V • $)r^ - (i^eL, 4>)r^ 

- (eu, Vi/zjTft 

- (e„, V ■ (0 0 tt)r^ + (i^eL, V</>)r^ - (ep, V ■ (i))r^ - (e„, y)r^ 

= - 1^(6^, V • na^)rh ~ (^cl, irG<i))7-^ - z/(e„, V • d$)rfe 

(Sixj ^^<3'0)7fc 

- (e„, V ■ (0 (g) + (Z^CL, V</>)r^ - (ep, V ■ z/zjr^ - (e„, 

taking {G,v,q,fi) = (z/I/g4), 0, I/q^, 0) in the error equation (I5.3I1 . inserting the resulting identity into the 
above expression and simplifying, we have 

= - (es, vnG^n)dTh - ^^(e«, V • 

- {(^u, nqilin)Qrn - (e^, Vd^)rfe 

- (e„, V ■ (0 (g ti)r^ + (i^eL, V</>)r^ - (ep, V ■ cl))r^ - (e„, y)r^ 

inserting two zero terms: (e^, v^rijQj-^ = (e^, ipn)Qjp and integrating by parts in the first two lines to obtain 

— (e^i e^, Ti -)- StpT ijQpp 

- (e„, V • {(t><Siu)r^ + (z/CL, V</>)r^ - (ep, V • (t>)r^ - {en,Y)r^. 

Next we work on the last line in the above expression. We first insert the projection of (p to have 

-(e„, V -{cptE) u)r^ + (z/eL, - (ep, V • - (e„, Y)r^ 

= iyei,,VIlvfpWh - (Spi ^ ■ nv(p)Th - (e«> V • (0 g u)r^ - (e„, Y)r^ 

+ (z/cl, '^S<i,)n - (ep, V • S^)r^ 
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taking {G,v,q,fi) = (0, 0, Z?m 0) in the error equation (|5.3L intergrating by parts for the last two 

terms in the above expression and simplifying, we have, 

= - e^), nv(f> - nM(f^)dTh + - 5pn — -iimSu, nv4> — 

- o{{u, It); (it. It), {nv(t),nM4>)) + o{{uh,Uh)\{uh,uh), {nv(t),nM4>)) 

- (e„, V ■ (</> (g) u)r^ - (e„, Y)r^ 

+ (UCL, ^S^)rh - (cp, V • S^)rh 

= —{-^{nMGu - Gu), nv(t> — nM(t>)dTh + - 5pn — -iimSu, nv4> — nM<t>)drh 

- ((e„, V • (</> (g) u)rh + e)((it. It); (e„, es), (-/7y^, HmcI)))^ 

- 0((ii, It); ((5„, 6u), {nv<f), nM(t>)) 

+ 0{{uh, Uh)] {Ufi, Uh), {nv(t>, nM4>)) - C>{{u, u); {Uh, Uh), {nv(f), nM(l>)) - (e„, Y)r„ 

+ {ve-L,^5,j,)Th - • ^ckWh- 

We can obtain the expression in the Lemma by inserting (^, </>) in the two O terms in the above identity. 
This completes the proof. □ 

Now we are ready to prove our last result: 

Lemma 5.5. Under the same assumption as in Lemma \5.‘A in addition we assume the full H'^ —regularity 
of the adjoint problem (lOl) holds and k > I, then we have 

lle.lln 

Here the constant C depends on ||it||fe+ 2 ,n, ||i^||ivuoo(n); ||p||fc+i,n,u and k but independent of h. 

Proof. By identity in Lemma 15.41 it suffice to estimate Ti — Tg. 

For Ti, we apply Cauchy-Schwarz inequality. Lemma 15^ (IS.lcI) and the regularity inequality (12.31) to have 

Ti < h~^ \\eu - eiaWdThl^^ \WS,s,n + S^nWoTh < Ch^+^ ■ /iH^(||$||yn + ||((<||i,n) < C/i''+^||e„||n. 
Similarly, for T 2 we have 

T2<i^h---\\nMeu- e^\\anh-"^\\nv<l>- nMHdn<Ch^^^ - h-'^h^Uh^n^Ch^^^WeJ^n. 

Using Cauchy-Schwarz inequality, (I5.1c|) . (I5.1dl) and (12.3L we can bound Tg as 

Tg < Ch^^^{\\L\\k+i.Q. + ||p||fe+i,n + ||it||fe+2,n) ll‘/’l|2,n < Ch^'^'^\\eu\\n- 

For Tg, we simply apply the Cauchy-Schwarz inequality, (I5.1bp . Lemma |5. 21 Lemma 15.31 and the regularity 
inequality (IQ) to have 

Ts < (u||eL||n + l|ep||n)||< ^^(ulleLlIn + ||ep||o)|| 0 || 2 ,n < C/i'"+^||e„||n. 

For Tg, we explicitly write this term: 

Tg = —{5u g> W, V77y^)7-h + {Tc{u){nMU — Hyu), Hvcf — nM(f)dTh + S^)n, nv4> — nM(f>)dTh 

< \\u\\oa,n(\\Su\\Th\\'^ nvffWn + \\nMU- f7yll||arfell^M0 - nv4>\\dTh + \\Su\\dTh\\nM4> - nycflldTh^ 

< C\\u\\^^Q{h'^+^\\u\\k+ 2 ,n\\(t>\\i,n + h'^~''^\u\\k+i,Qh^\4>\\2,n) by (I5.1b|) . (15.1cD and (I5.1dl) 
<C/i'=+2||e„||n, 


by the regularity assumption (EH). 
















16 


WEIFENG QIU AND KE SHI 


For T4, we first expand the term as: 

Ti = - (e„, (V • - (cu ® u, V4 >)t^ + (e„ 0 u, vnY(f))r^ 

— {Tc{u){eu — e^), nv4> — nM(l>)dTh ~ ii^u ® u)n, nv4> — nM(f>)dTh 
= — (e„ 0 M, — {Tc(u-)(eu — e^), nv4> — nM(l>)dTh ~ (i^u ® Uycf) — nM(t>)dTh 

<C'||«||oo,n(||e„||n||Vd0||rh +h~^\\eu- - nM(t>\\dTh + WfiuWdThW^vffi - nM(t>\\dTh) 

<C/i'"’+^||e„||n + Ch^ ||e„||n||e^||arhj 

by Lemma [5TH (jS.lbl) . (j5.1cp and (IS.ldl) . By a triangle inequality we have 

WduWdTh < ll^ix “ (iuWdTh + l|e»||arft < ). 

Inserting this inequality into the estimate for r4 we obtain: 

l4<C/i'=+2||e„||n. 

To bound Tq, we first derive some useful inequalities, we first bound ||ith||oo,n: 

(h.8) ||i*/i||oo,si ^ ll^ulloo,!! “f ||-ffv/ii||oo,o ^ (J{h 2 ||e^||n T ||'ii'|joo,o)- 

Next by a triangle inequality, we have 

(5.9) \\uh — u/illaTfi ^ W^u ~ ^uWdTh + W^vu— IlMuWdTh ^ . 

Consequently, we have 

(5.10) \\Uh\\oc,dTh ^ Wu-h — UhWao^dTh + \\''J'h\\ao,dTh < ^ + C {h~ ^\\eu\\n + \\u\\ao,Q) ■ 

The last step we applied a scaling argument for the polynomials on dTh- Finally, applying triangle inequality 
we obtain the following estimates: 

(5.11a) ||tt- Uh\\n < ||e„||n + HiJkIIo < 

(5.11b) ||V(1X- t.,)||r, < ||Ve„||r. + ||Vd„||r. < 

(5.11c) \\u-Uh\\dTh. < ll^sllaTh + lle^llar^ 

Now we are ready to present the estimate for Tg. If we expand Tq using the definition of O, we obtain: 

Te = - {uh 0 {u- Uh),yS^) - (i(V • {u- Uh))uh,S^) + {]^{uh 0 (w/i - Uh))n, 5(f,)dTh 
+ {{tc{u) -Tc{uh)){uh - Uh),6^ - 5-^)dTh + {{^h ® {u-Uh))n,6^ - S^)oTh 
applying the generalized Holder’s inequality for each term, we have 

^ II ||oo,Q (11^- Wh||a||V(5<^||r;, + ||V • {u- Uh)\\rh\\^(i>\\n) + HwhUoo.nllw/i - Uh\\dn\\^(i>\\dTh 

+ (l|M||oo,a + \\Uh\\oo,dTh)\\'>J-h - Uh\\dTh\\S<l> - + \\Uh\\oc,dTh\\u - Uh\\dTh¥(l, “ S^\\oTh 

now if we apply the inequalities (15.81) - (I5.11|) . (IS.lbl) - (IS.ldl) and (12.31) . we have 
<C/i'=+2(h'=+i-i + l)||e„||n < C/i'=+2||e„||n. 

Finally, we need to estimate T7 which is more involved than the previous terms. To this end, we begin by 
expanding the nonlinear operator O: 

Tr = {uh 0 {u- Uh),V4))T^ + (iv • {u- Uh)uh,4>)r^ + {^{uh 0 (tt/i - Uh))n,(f>)or^ - {eu,Y)r^, 
integrating by parts the second term, we have 

= {uh 0 (m- Uh),Vcf))r^ + {^{Uh 0 (u- Uh))n, 4>)ort, - i^Uh iSi{u- Uh),V4>)r^ 

- {u-Uh),Vuh)Th + {^{Uh 0 {uh-Uk.))n,(f))on - ieu,Y)r^ 

= 0 {u- Uh),Vcf))rf^ + 0 {u-Uh))n, 4>)drH - {i^4>®{u-Uh),Vuh)Tt, - (e„,y)r;,. 














A SUPERCONVERGENT HDG FOR THE INCOMPRESSIBLE NAVIER-STOKES EQUATIONS 


17 


inserting the zero term —{u — Ui^))n, 4 >)dTh — 0 above expression, 

= Uh) ® {u-Uh))Ti,(f))on + ® {u- Uh),V(f>)r^ - {u-Uh),Vuh)rH. - (e«,l")rh 

= Uh) ® {u-Uh))Tl,(f))on - {^{u-Uh) ® {u- Uh),V4>)r^ + {u- Uh),y{u- Uh))Th 

+ {u- Uh),V(t>)n - {^4>® (u- Uh),Vu)rh - {eu,Y)rh, 

by the definition of Y = we obtain: 

= Uh) ® {u-Uh))Ti,(f))ou - {^{u-Uh) 0 {u- Uh),\/4>)r^ + {u-Uh),y{u-Uh))Th 

+ {^u0Su,V(t>)r^ - {^(f>®dy,,Vu)rh 

= T 71 + ■ • ■ + T 75 . 

We are going to estimate each of the above terms. For T 71 we apply the generalized Holder’s inequality, 
(I3.1c|) . (15.111) . and (12.3L 

T71 < ||tt- Uh\\L‘^{drh)\\'^~^h\\dTh\\(t^\\L'^{dTh) ^ Ch~^\\u- Uh\\l,h\\<t>\\l,h\\'^ ~ ^hWdTh < ||eu||n. 

For T 72 , we apply the generalized Holder’s inequality, p.lap . (15.lip and (12.31) to get: 

T 72 < ||m- u/.||i4(o)||V0||n < c\\u- 

Similarly, we can bound T 73 as 

T 73 < ||^||L'‘(a)||M- Ufi\\L^(Q)\\V{u- Uh)\\Th < C/i^''+^||e„||a. 

For T 74 ,r 75 we apply the generalized Holder’s inequality as 

T74 < ||«||oc.n||^K||n||V<A||a < C||M||oo.a/i'=+'||e„||a, 

T 75 < ||Vu||oo.n||<A||o||<5„||o < C||V'a||oo.n/i'=+'||e„||n. 

The proof is complete by combining all the estimates for Ti — Tg. □ 

6 . Concluding remarks 

In this paper, we introduced and analyzed a new HDG method for the Navier-Stokes equations. The 
work can be seen as a continuation of our previous work on HDG methods for linear problems, see [25l 126) . 
Gomparing with the original HDG method for Navier-Stokes equation ISlEs], our method uses an enriched 
polynomial space for the velocity in each element, a modified numerical flux and a modified HDG formulation. 
As a concequence, we obtained optimal order of convergence for all unknowns and superconvergence for the 
velocity without postprocessing. In addition, similar as in [251126] , the analysis in this paper is valid for 
general polygonal meshes. 

Numerical study of the method as well as other computational aspects including the characterization of 
the scheme, implementation of the method using Picard iteration and other related issues will be extensively 
discussed in a separate paper. 
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